Ab-initio No-Core Gamow Shell Model calculations with realistic interactions 
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No-Core Gamow Shell Model (NCGSM) is applied for the first time to study selected well-bound 
and unbound states of helium isotopes. This model is formulated in the rigged Hilbert space and, by 
using a complete Berggren ensemble, treats bound, resonant, and scattering states on equal footing. 
We use the Density Matrix Renormalization Group method to solve the many-body Schrodinger 
equation. To test the validity of our approach, we benchmarked the NCGSM results against Faddeev 
and Faddeev- Yakubovsky exact calculations for 3 H and 4 He nuclei. We also performed ab initio 
NCGSM calculations for the unstable nucleus 5 He and determined the ground state energy and 
decay width, starting from a realistic N 3 LO chiral interaction. 
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I. INTRODUCTION 

In the last decade our knowledge of nuclei far from 
the valley of stability has radically improved. This im- 
provement is a by-product of advances in both experi- 
ment and theory. New experimental facilities that have 
already been built (RIBF at RIKEN) or are being con- 
structed (SPIRAL2 at GANIL, FAIR, FRIB at MSU) will 
give us a better insight of areas in the nuclear chart that 
have never been explored, pushing even farther the lim- 
its of nuclear existence. A few decades ago, the nuclear 
chart consisted of approximately 1000 isotopes, whereas 
in 2011 this number has been expanded to approximately 
3000 species, and an estimated number of nuclei that can 
exist in nature or synthesized in the laboratory is approx- 
imately 7000 [l|. The increase in computing power has 
made it possible to calculate the properties of nuclei in 
an ab initio manner, using realistic interactions, which 
reproduce the nucleon-nucleon scattering data. For few- 
body systems (A < 4) methods such as Faddeev Q and 
Faddeev- Yakubovsky (FY) [3j provide an exact solution 
to the many-body problem. Methods such as the Green's 
Function Monte Carlo (GFMC) [3], the Hyperspherical 
Harmonics Hi], the No-Core Shell Model (NCSM) 
the Coupled-Cluster approach (CC) @, and more re- 
cently, the In-Medium Similarity Renormalization Group 
method (IM-SRG) [T3, EH and Dyson Self- Consistent 
Green's Function (Dyson-SCGF) method [HI have been 
applied successfully for the ab initio description of light 
and medium mass nuclei. 

Nuclei with a large isospin that can be found in re- 
gions far away from the valley of stability have attracted 
a great deal interest. They belong to the category of 
Open Quantum Systems (OQS) 13], which in the case 
of the nucleus are inter-connected via the decay and 
reaction channels. They are very fragile objects with 
small separation energies and very large spatial dimen- 
sions. The proximity of the continuum affects their bulk 
properties (matter and charge distributions) and their 



spectra. Phenomena such as the anomalous behavior of 
elastic cross-sections and the associated overlap integral 
near threshold states in multi-channel coupling (Wigner- 
cusps) [H],[lij], the isospin and mirror symmetry-brea king 
threshold effects [HI, [uj, the resonance trapping []~8p - 
|20( and super-radiance phenomenon [2l|, |22j , the appear- 
ance of cluster correlations in the vicinity of the respec- 
tive cluster emission threshold [23|, the modification of 
spectral fluctuations [24T - |26| , and deviations from Porter- 
Thomas resonance width distribution [2(| [27], HH , are all 
unique manifestations of the continuum coupling. 

For their theoretical explanation it was necessary to 
generalize existing many-body methods, and create the- 
ories which unify structure and reactions. Examples of 
these attempts are the Shell Model Embedded in the 
ContinuumisMEC) [H-[3l| and the Gamow Shell Model 
(GSM) [H-il- The SMEC is a recent realization of the 
real-energy Continuum Shell Mode l [3a . l36j which uses 
the Feshbach projection technique [37j in order to take 
into account the coupling to the scattering continuum. 
The GSM is a generalization of the Harmonic Oscillator 
based shell model in the complex energy plane by using 
the Berggren ensemble [38], which treats resonant and 
non-resonant states on equal footing. Ab initio calcu- 
lations that use a proper asymptotic behavior of wave 
functions include the NCSM coupled with the Resonat- 
ing Group Method [H, Ho| and the CC approach gen- 
eralized in the complex-energy plane using the Berggren 
basis [Mlil- 
In this work we introduce the No-Core Gamow Shell 
Model (NCGSM) as an alternative for calculations of 
weakly bound and unbound states of light nuclei using 
realistic interactions and allowing all the nucleons to be 
active. The paper is organized in the following manner: 
in Sections |TT] and IIII1 we describe the basic ingredients 
of our method, such as the many-body Hamiltonian, the 
single-particle basis we employ, the way the two-body 
matrix elements are calculated within the Berggren ba- 
sis, and we discuss the translational invariance of our 
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Hamiltonian. In Section HVl we describe the Density Ma- 
trix Renormalization Group (DMRG) method, which is 
an efficient tool for a diagonalization of large complex- 
symmetric GSM matrices. In Section [V] we present our 
calculations for the 3 H, 4 He and 5 He nuclei and, finally, 
in Section fVTl we discuss the conclusions and the future 
perspectives. 



II. HAMILTONIAN 

Our goal is to solve the ^4-body Schrodinger equation 
H\*) = E\9), (1) 
where H is the intrinsic Hamiltonian 
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and V NN a realistic NN interaction. For ([2]) the following 
identities are useful: 
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FIG. 1. (Color online) An illustration of the Berggren s.p. 
basis used in the NCGSM, showing the position of resonant 
(bound states and resonances) states in the complex fc-plane. 
The non-resonant continuum states lie along the complex con- 
tour L+. 
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where P is the Center of Mass (CoM) momentum, and 
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In the NCGSM there is no restriction on the type of the 
NN interaction, contrary to the GFMC approach for ex- 
ample, where difficulties arise for non-local potentials. 
One can use a local interaction, such as the Argonne v\% 
potential liH. o r a non-local interaction, such as the CD- 
Bonn 2000 [44( and various chiral interactions. There is 
also the possibility to use renormalized versions of the 
aforementioned forces, by applying techniques such as 
Viowk [451 ] . the Similarity Renormalization Group (SRG) 
approach [H, H3], or the G-matrix [48]. In this work 
we employ the phenomenological Argonne Wig potential 
and the chiral N 3 LO interaction [49] which is consistent 
with the symmetries of the QCD Lagrangian. Both po- 
tentials were renormalized via the Vi ow & method with a 
sharp momentum cutoff A = 1.9 fm -1 to decouple high 
from low momentum degrees of freedom and, henceforth, 
improving the convergence of nuclear structure calcula- 
tions Moreover, specific interactions will be used to 
compare NCGSM with other approaches. 



III. BERGGREN BASIS 

In previous applications of the GSM, where a tightly 
bound core was assumed ( 4 He or 16 0), the single parti- 
cle (s.p.) basis was usually generated by solving the one- 
body Schrodinger equation with a Woods-Saxon (WS) 
potential, parameterized to reproduce the core plus nu- 
cleon spectrum. In the case of the NCGSM, the s.p. 
basis will be generated by the realistic two-body interac- 
tion itself by solving the integro-differential Schrodinger 
equation which contains both local and non-local parts 
[HH HH . This numerical method is known as the Gamow 
Hartree-Fock (GHF) method since it generates a mi- 
croscopic basis that includes resonant and non-resonant 
states. The GHF method can be applied not only in 
spherical cases (closed-shells) but also in deformed cases 
(non-closed shells) 51]. The one-body self-consistent 
potential \Jhf(y) is then used to solve the one-body 
Schrodinger equation: 
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— + ^(UHF(r)+V c (Z c ,r))-k 2 



where V c (Z c ,r) is the one-body Coulomb potential: 
V%Z c ,r)= CcZcerf{ar) 



Uk(r) 
(6) 

(7) 



and C c is the Coulomb constant, Z c the proton number 
and a is a constant, which is given by a = 3^/tt/ARq. 
The reason we choose an error function to approximate 
the Coulomb field and not for example, the field that is 
produced by a uniformly charged sphere at Rq, lies in 
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the fact that the latter is non-analytic at Rq. The value 
of Rq is chosen in a way that the Coulomb potential of 
Eq. (|7|) and the potential of a uniformly charged sphere 
are equal at the origin. The wave number k is defined 
as: k — \j2mE /h and is, in general, complex. Equation 
^ is solved with the requirement that at large distances 
the wave function will behave as a linear combination of 
Hankel or Coulomb functions for neutrons and protons, 
respectively: 

u fc (r) ~ C + H ( e +\kr) + C~ H$ '(kr) (8) 

where r) is the Sommcrfcld parameter. The C + and C~ 
coefficients are determined by the normalization of the 
radial functions u(r) to a Dirac's 5 distribution: 

POO 

/ u k (r)u k ,(r)dr = S(k-k') (9) 
Jo 

The solutions of ([5]) which satisfy pure outgoing bound- 




r(fm) 



FIG. 2. (Color online) Selected squared radial basis func- 
tions in the NCGSM. Resonant (bound, resonances) and non- 
resonant (scattering) states are plotted. The basis generating 
potential is produced by the Vi ow k N 3 LO interaction for the 
5 He nucleus. The notation nlj is explained in the text. 

ary conditions (C~ — in ©) correspond to the poles 
of the S-matrix and they are represented as dots in the 
complex fc-plane of Fig. [1] While normalization of bound 
states does not pose any difficulty, one should pay more 
attention on the normalization of resonances. The latter 
diverge exponentially for large distances and the regular- 
ization method that is used for the calculation of their 
norm is the external complex scaling 34]. The method 
facilitates from the fact of using complex radii for the 
integration of resonant wave functions: 

f oo pR poo 

/ u 2 k {r)dr= / u\{r)dr + / u k {R + xe l9 ) 2 e w dx, 

J0 Jo J R 

(10) 

with R chosen sufficiently large so as to match the con- 
dition ([SJ for C~ = 0, while 9 is the angle of external 



rotation, which satisfies the condition that u k (R + xe l9 ) 
= for x — » oo. 

It was shown by Berggren [38| that for a given partial 
wave (£,j) the scattering states, which are distributed 
along the L + contour and the resonant solutions (bound 
states and/or resonances) of © form a complete set: 

^\Un){u n \+ \u k ){u k \dk=l. (11) 

n ^ L + 

The completeness reassures us that any function which 
lies between the contour and real fc-axis and exhibits 
outgoing wave asymptotic (e ikr ), can be expanded us- 
ing (TTTI) . In practice the integral in (fTTj) is discretized 
by means of an appropriate quadrature rule (the Gauss- 
Legendre quadrature in our case) and we end up with a 
discretized completeness relation: 

N 

2jw n |ix„)(u„| ~ 1, (12) 

i=l 

where uj n — 1 both for resonant states and the Gauss- 
Legendre weight for non-resonant states along the dis- 
cretized contour. The approximate equality in (fT2|) arises 
from the finite discretization of the contour. In addition, 
the discretized contour does not extent to infinity and 
we use a maximum cutoff (k max ) at around 4 fm" 1 . We 
have checked that results do not depend on the values of 
the k m ax for an adequately large number of discretization 
points. 

In Fig. [2]we show the radial behavior of a few Berggren 
basis states that were generated from the N 3 LO interac- 
tion with a Viowk cut-off A = 1.9 frn -1 . The states that 
are plotted refer to neutron states. The 0^3/2 resonant 
state has a complex energy and lies in the fourth quad- 
rant of the complex fc-plane (see also Fig. [T|). It is a solu- 
tion of Eq. ^) with outgoing wave boundary conditions 
at large distances. We, indeed, observe a localization of 
the state in the region of the attractive nuclear potential 
(r < 2 fm) and an outgoing wave behavior for distances 
beyond the range of the nuclear potential. This is the 
basic characteristic of a metastable s.p. state. 

In order to satisfy Berggren's completeness, the scat- 
tering contour L + has to be complex. It is then un- 
derstood that the 10p 3 / 2 is a point along the complex 
contour. It corresponds to a state, which is given as a 
linear combination of Hankel functions, as it is seen in 
Eq. ([5]). Here we are using the notation n£j, where n is 
the radial quantum number and is identified as the 10th 
Gauss-Legendre discretization point on the L + contour; 
I is the s.p. angular momentum of the state (£ = 1 in 
this example), and j is the s.p. total angular momen- 
tum 0. The 0si/2 resonant state is bound and lies on 



For simplicity we assume partial waves with one pole. The pole is 
always identified with the n=0 shell in our algorithm and, hence, 
the scattering states start from the (n = 1) shell. 



4 



the imaginary momentum axis with a real and negative 
energy. At large distances, the state state is decaying as 
an exponential. Also shown are the 10si/2 and 10pi/2 
scattering states, which both lie on the real momentum 
axis, since there is no complex resonant state associated 
with these partial waves. Similar to the complex 
state, they correspond to the 10th discretization point 
of the Gauss-Legendre rule, but the contour lies along 
the real momentum axis. At this point we would like to 
highlight that the Berggren basis not only imposes the 
correct quantum mechanical asymptotic behavior of the 
s.p. states but also includes the continuum states in a 
rigorous manner, promoting it to an ideal and realistic 
basis for the description of metastable and weakly bound 
states. 

The completeness relation (TT2"j) is similar to a usual 
discrete completeness relation, such as the Harmonic Os- 
cillator (HO)) one, and results in an eigenvalue prob- 
lem. The many-body basis states are Slater determinants 
\SD n ) = \u\ . . . ua), where Uk is a resonant (bound state 
or resonance) or non-resonant (scattering) state. In this 
basis, the Hamiltonian matrix is complex symmetric and 
upon diagonalization, many-body correlations and cou- 
pling to the continuum are taken into account simulta- 
neously. As a direct consequence of Eq. (fT2|) . the many- 
body states also satisfy the completeness relation: 

N 

J2\SD n )(SD n \~l. (13) 

i=l 

The squares of the expansion coefficients, which are de- 
termined by the diagonalization and not by their absolute 
values, satisfy the relation: 

JV 

£4 = 1. (14) 

i=l 

Furthermore, the completeness relation (I13|) can be 
used for the calculation of two-body matrix elements 
(TBMEs) between Berggren basis states. 



A. TBMEs of realistic interactions in the Berggren 

basis 

— Nuclear part — Matrix elements of realistic interac- 
tions are defined in a relative and CoM system of coor- 
dinates. In order to work in a basis of Slater determi- 
nants, a transformation from the relative and CoM to 
the laboratory system is necessary. When working in 
the HO s.p. basis, this is possible through the Brody- 
Moshinsky brackets (53|. For a different basis such as 



2 This can be immediately checked by the fact that a bound reso- 
nant state (e lkr ) behaves as e — Kr for r — > oo, since k = in. 



the Berggren basis, one has to perform a multiple de- 
composition of the realistic NN interaction and calculate 
two-dimensional radial integrals. Matrix elements be- 
tween scattering states need to be regularized by means 
of the complex scaling (CS) technique [34| which unfor- 
tunately does not work for just any type of integrals and 
can cause numerical instabilities. The problem is allevi- 
ated by expanding the NN interaction in a truncated HO 
basis |54| : 

V NN = \ap)(ap\V NN \75)(^5\ (15) 

The matrix elements of the NN interaction in the 
Berggren ensemble become then: 

(ab\V NN \cd) = £ (ab\aP)(a(3\V NN \~f6)(~f6\cd) (16) 

We end up calculating overlaps between HO and 
Berggren states (af3\ab) and (jSlcd), where the Latin let- 
ters denote Berggren states and Greek letters HO states. 
Due to the Gaussian fall-off of the HO states, no complex 
scaling is needed for the calculation of these integrals. 
On the other hand, matrix elements of the NN interac- 
tion in the HO basis (aj3\VNN\l5) can be conveniently 
calculated using the Brody-Moshinsky brackets [53[. 

The method of handling matrix elements in the 
Berggren basis using the projection of continuum onto 
the HO states should not be confused with basis expan- 
sion methods suitable for the description of closed quan- 
tum systems. Only the short range part of the nuclear 
interaction is expanded in the HO basis. The kinetic en- 
ergy operator is calculated in the Berggren basis, so the 
calculations of weakly bound and unbound systems are 
possible. With this formulation it is also clear that there 
is no restriction on the type of NN interaction one can 
use in the NCGSM. What we need is just the nuclear 
matrix elements in the HO basis. 

— Coulomb interaction treatment — For nuclear sys- 
tems with two or more protons the two-body Coulomb in- 
teraction is included in the Hamiltonian @ . The method 
we adopt for the treatment of the long-range Coulomb 
interaction was first used in the description of isospin 
breaking due the continuum coupling [17| and it was also 
recently applied to calculate reaction observables [H, H|[ 
(see also [53] for a detailed numerical analysis). The ba- 
sic idea of the method is to add and subtract the one- 
body Coulomb potential ([7]) (with Z = 2, e.g., in case 
of 3 He or 5 He) from the two-body Coulomb interaction: 
K(l,2) = U c (l) + 04(1,2) - £4(1)). Then the second 
term in the parentheses has a short-range character and 
the HO expansion method of (|T~5]) can be applied. Matrix 
elements of the Coulomb interaction can be calculated 
using the Brody-Moshinsky brackets without the need to 
perform an external CS calculation. 

We would like to mention here that this method of 
treating the two-body Coulomb interaction would be of 
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particular interest, when one has to deal with many-body 
proton resonances, such as in the 6 Be nucleus. In this 
case, calculating the Coulomb interaction potential by 
simply expanding it in a HO basis (see Eq. (1151) ). would 
be a rather poor approximation. 

— Center of mass ( CoM) motion in the NCGSM — By 
adopting a s.p. basis upon which we build many-body 
basis states, we effectively localize the nucleus in space 
and, hence, we break the translational invariance of the 
Hamiltonian. Moreover, we would need ZA — 3 coordi- 
nates to describe it, but the nuclear wave function we 
construct, depends on 3^4 coordinates, where A is the to- 
tal number of particles. These redundant degrees of free- 
dom are responsible for the CoM spuriosity that appears 
in many-body methods. On the other hand, plane waves 
s.p. states are eigenstates of the momentum operator 
and, hence, preserve the transitional invariance, but un- 
fortunately cannot be used to describe a localized system. 
The alternatives are: i) Solving the many-body problem 
using relative coordinates (e.g. Jacobi coordinates) which 
reassures the translational invariance of the system, with 
the price of unfeasible antisymmetrization of states for 
A >8, and ii) Using the unique analytical properties of 
the HO s.p. basis in a full Nfku space, in which the total 
wave function is factorized into |'0rei)<8>|?/'CoM), limiting 
though the application to well-bound systems only. In 
the case of the NCGSM the latter factorization is not 
guaranteed and it has to be demonstrated numerically. 
Since our Hamiltonian @ is intrinsic, we are expecting 
that in an infinite space there would be no spuriosity. 
However, because we are working in a finite space, it is 
necessary to check numerically this condition. 

Assuming that the factorization into a CoM and a rel- 
ative wave function is valid and also the CoM wave func- 
tion has a Gaussian shape, we calculate the expectation 
value of the CoM operator [58| : 

1 -j 2 mAu 2 -> 2 3 
HcoM = ~^A P cm H 2 — cm ~ 2 ' ' 

where hw is the parameter that characterizes the Gaus- 
sian wave function. The matrix elements of (|17[) are cal- 
culated with the HO expansion method of (fTS)) and the 
analytical formulas for their expression are found in [59j j . 
Following the assertion of Rcf. 60], if (H cm ) ~ then 
the factorization is valid. 



IV. RESOLUTION OF THE MANY-BODY 
SCHRODINGER EQUATION WITH THE DMRG 
METHOD 

The Schrodinger equation (|TJ) is solved within a many- 
body basis constructed from the discrete set of single- 
particle states \ui) (fT2]) . The discretization of the inte- 
gral along the contour L + in Eq. (|TT|) should be precise 
enough so that the discretized completeness relation (fT"2")) 
is fulfilled. In other words, the number of discretized 



shells should be increased until Eq. (|12[) holds. As a con- 
sequence, the dimension of the many-body model space 
will increase dramatically with the number of nucleons 
and number of shells. Efficient numerical methods al- 
lowing the diagonalization of large Hermitian as well as 
complex-symmetric matrices are then required to solve 
the NCGSM problem. 

In this paper, we have used one of these methods, 
namely, the DMRG method [6l|, [62| which has been gen- 
eralized in the context of the GSM in [H 11]. The 
GSM/DMRG approach has been applied previously to 
study several weakly-bound/unbound nuclei described 
as few-valence-nucleon systems interacting via schematic 
two-body forces above an inert core. In this paper, all 
nucleons are considered active and realistic two-body in- 
teractions are used but nevertheless, the application of 
the DMRG method is similar. In the following, we re- 
call the main ideas of the DMRG in the multi-shell GSM 
problem |64j |. 

The purpose of the DMRG method is to allow the cal- 
culation of the many-body poles of the scattering ma- 
trix of the NCGSM Hamiltonian H by performing ef- 
ficient truncations of the many-body model space. As 
the contribution of the non-resonant continuum to the 
structure of many-body bound/resonant eigenstates of 
H is usually smaller than the contribution from the 
bound/resonant orbits, the following separation is per- 
formed: the many-body states constructed from the s.p. 
poles form a subspace H (the so-called 'reference sub- 
space'), and the remaining states containing contribu- 
tions from non-resonant shells form a complement sub- 
space P. The set 8 of many-body basis states (fT3"|) can 
then be written as 

£ = H®P. (18) 

The DMRG technique is then used to performed trunca- 
tions in £ by keeping only selected 'optimized' states in 
P in the sense of a criteria based on the density matrix 
in P (see below for a more explicit explanation). 

Let us assume that we want to calculate an eigenstate 
l^) of H for a nucleus coupled to the total angular mo- 
mentum J and parity 7r. The number of proton(s) and 
neutron(s) are respectively and n„. One begins by 
constructing all states \k)n forming the subspace H . The 
set of those states is denoted as {£;#}. The many-body 
configurations in H can be classified in different families 
{n; jfj} according to their number of nucleons n, total an- 
gular momentum jn, and parity ir. States with a number 
of protons (neutrons) larger than (n„) are not consid- 
ered since they do not contribute to the many-body states 
in the composition of subspaces H and P. The matrix 
elements in H of the suboperators of the NCGSM Hamil- 
tonian H expressed in the second quantization form, are 
calculated and stored: 

{0} = {a\ (a f a) K , (aV) K ((a t a t ) A '5) L , 

(a^a)) K {aa) K }, (19) 
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where a* and a are the nucleon creation and annihilation 
operators in shells forming the subspace H . The NCGSM 
Hamiltonian is then diagonalized in the pole space H to 
provide the zeroth-order approximation |\&}( ) to |^). 

In the following stage, the subspace P is built, step 
by step, by adding scattering shells one by one during 
the so-called 'warm-up phase'. At each step, many-body 
states constructed within the new added shell are cou- 
pled to optimized many-body states constructed during 
previous steps i.e., constructed within previously added 
shells. Moreover, the matrix elements (fT§|) of the sub- 
operators acting among the optimized states have been 
stored during previous steps. 

To be more specific, let us assume that the s th step 
is reached. The method is illustrated in Fig. [3] The 
scattering shell (lj) s belonging to the discretized con- 
tour L + is added and within this shell, one constructs 
all possible many-body states {{lj)^ p }- Matrix elements 
of suboperators (IT!?)) acting on {{lj)™ p } are also com- 
puted. One then couples previously optimized states de- 
noted as \ot)p to {{ljYs F } to obtain the set of states 
{ip} = {ap ® Gi)s P }- States in H are then coupled 
with the states \ip) to construct the set {kn ® ip} J of 
states coupled to J 71 ' that constitutes a basis in which the 
NCGSM Hamiltonian is diagonalized. The Hamiltonian 
matrix is constructed in this set with the Wigner-Eckart 
theorem and the matrix elements of the suboperators 
(OH) acting on {k H },{a P } and {(/j)" p }- 

The target state is selected among the eigenstates 
of H as the one having the largest overlap with the ref- 
erence vector l^)' ^. Based on the expansion 
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EKhUh 
ip 

k H ,i P 



'l${\k H (j H ))®\i P (jp))} J , (20) 



by summing over the reference subspace H for a fixed 
value of jp , one defines the reduced density matrix (65j : 



Pi 



'ip(Ji 



(21) 



Truncation in the subspace P is dictated by the den- 
sity matrix. In the case of a Hcrmitian Hamiltonian, the 
eigenvalues of the density matrix are real and one can 
show that the truncation in P is optimal when one keeps 
the eigenstates of the density matrix with the largest 
eigenvalues j§6| . More specifically, the error in the repre- 
sentation of \ip) (I2TJ1) after truncation, is minimal in that 
case. In that sense, the eigenstates of the density matrix 
are 'optimal'. 

Within the metric defining the Berggren ensemble, the 
NCGSM density matrix is complex-symmetric and its 
eigenvalues are, in general, complex. The trace of the 
density matrix being equal to one, the truncation is done 
by keeping eigenstates of the density matrix with the cor- 
responding eigenvalue w a such that the condition 




< e 



(22) 



{(lj)s} 

{a P } 



FIG. 3. Schematic illustration of the NCGSM/DMRG pro- 
cedure during the s th step of the warm-up phase. States 
{kn} from H, previously optimized states cxh, and states 
{{lj) s } constructed by occupying the s th shell with n parti- 
cles are coupled to generate the new set of states {kn®ip} J = 
{k H ® {a P ® (Zj)"}} J - 



is satisfied. The quantity e in (f2"2"|) can be viewed as 
the truncation error of the reduced density matrix. The 
smaller e, the larger number of eigenvectors must be 
kept. In particular, for e=0, all eigenvectors with non- 
zero eigenvalues are retained. 

One then keeps eigenstates of the density matrix ac- 
cording to Eq. (f2"2"|) . These are expressed as linear com- 
bination of the vectors \i)p in P and all matrix elements 
of the suboperators in these optimized states are recal- 
culated and stored. Note that at each step, we enforce 
that at least one state in each family {n;jp} is kept [64j . 

The warm-up phase continues by having the P sub- 
space grow by adding scattering shells one by one until 
the last shell is reached, providing a first guess for the 
wave function of the system in the whole ensemble of 
shells. At this point, all s.p. states have been consid- 
ered, and all suboperators of the Hamiltonian H acting 
on states saved after truncation in P have been computed 
and stored. The warm-up phase ends and the so-called 
sweeping phase begins. 

Starting from the last scattering shell (Ij)iast, the pro- 
cedure continues in the reverse direction (the 'sweep- 
down' phase) using the previously stored information. 
At this stage, the truncations are done according to the 
density matrix, which is obtained by summing over states 
\kn) of the reference subspace H and the states \i pre v) 
generated in the warm-up phase. Scattering shells are 
added one at a time and at the last step of the sweep 
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down phase, the first scattering shell is reached. The 
procedure is then reversed and a sweep in the upward 
direction (the 'sweep-up' phase) begins. Using the infor- 
mation previously stored, a first shell is added, then a 
second one, etc. The sweeping sequences continue until 
convergence for the target eigenvalue is achieved. For 
more details, see Ref. |64| . 



gence and in all the following we will adopt n„ 
for our calculations. 
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B. Convergence with respect to the angular 
momentum of the model space and with respect to 
the number of sweeps in the DMRG method 



V. APPLICATIONS OF THE NCGSM/DMRG 

A. Test of convergence with respect to N max of the 
two-body interaction in HO expansion 

In Table U we check the energy convergence of the 3 He 
nucleus with respect of the N max parameter in the ex- 
pansion of the two-body interaction. We choose 3 He due 
to the existence of both nuclear and Coulomb parts in 
the NCGSM Hamiltonian. The interaction employed is 
the N 3 LO renormalized at a cut-off A = 1.9 fin -1 , and 
we used a s.p. Berggren basis consisting only of s±/2, 
Pi/2 and p 3 /2 orbitals, for both neutrons and protons. 
The Osi/2 states are bound and the scattering contours 

TABLE I. Convergence of the NCGSM eigenvalues with re- 
spect to the N m ax = "2n ma x + imax of the HO expansion of the 
NN interaction. Values are in MeV. The length parameter of 
HO states is b = 1.5 fm. For this value of the b parameter, 
the CoM wave function is a Gaussian. 



max 


Energy 


5 


-5.321 


7 


-5.334 


9 


-5.336 


11 


-5.343 


13 


-5.343 


15 


-5.343 


17 


-5.346 


19 


-5.349 


21 


-5.352 


23 


-5.352 


25 


-5.352 


27 


-5.352 


29 


-5.352 



lie on the real momentum axis. They are discretized 
with 20 points each and they extend up to 4 fm -1 . At 
this point, our purpose is not to describe realistically 3 He 
but rather check the convergence of the many-body result 
with respect to a number of HO states of the expansion 
of the realistic interaction. In Table U we see that there 
is an overall small variation of the energy with increasing 
(or N max ). The energy changes only by ~ 30 keV 
by changing N max from N max = 5 to N max = 29. For 
N m ax = 21 {n max — 10), the energy shows the conver- 



In this section, we are applying the NCGSM method 
for the 3 H nucleus. 3 H (and also 4 He or 3 He) are well 
bound systems and the HO basis is sufficient for their 
description. Namely, the continuum has a negligible im- 
pact on their g.s. properties. It is then understood that 
our purpose is to solely test our approach against well- 
known bound-states methods. 

In Fig. 21 we compare NCGSM results obtained with 
the Argonne vi$ against Faddeev calculations. We per- 
form a calculation using s 1/2 , P3/2, P1/2, d 3/2 , 4/2, A/2, 
/5/2, 59/2, 57/2 partial waves for protons and neutrons. 
The basis generating potential is the GHF which gives 
the Osi/2 proton and neutron states bound, with ener- 
gies -10.417 MeV and -11.982 MeV along the imaginary 
momentum axis. For this reason, the scattering continua 
*{ s i/2}, i{P3/2},i{Pi/2} are chosen along the real axis and 
each of them is discretized with 18 points. Here i denotes 
the number of discretization points which ranges from 1 
to 18 for the s-wave and from to 18 for the p- waves. For 
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FIG. 4. (Color online) 3 H results as a function of the (maxi- 
mum) angular momentum, in the interval from I = to I = 4. 

the d- waves we choose a discrete HO basis with b = 1.5 
fm. In particular, we use five HO states for the d%/2 an d 
d§/2 states, which from now on we denote as 5d. The 
physical argument behind this choice lies in the fact for 
I > 1 the centrifugal barrier is large enough to confine the 
s.p. states and, hence, the HO basis becomes a realistic 
alternative. For / and g partial waves we used three HO 
states with b = 1.5 fm for both protons and neutrons. In 
total, the s.p. basis consisted of 154 partial waves. The 
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result we obtain is: 

£ncgsm = -8.39 MeV, 
whereas the Faddeev result [6tJ is: 

^Fadcov = -8.47 MeV. 

In FigHJ we show also an exponential fit of the NCGSM 
results for different £ of the s.p. basis. The extrapolated 
result is 



E, 



cxtrp 



= -8.449 ±0.087 MeV, 



where the fit function is: E = E, 



cxtrp 



b ■ e 



-c-l 



We see 



that inclusion of partial waves with angular momentum 
larger than £ = 4 has a very small contribution of the 
order of 50 keV to the many-body result. 

For three-particles systems we can perform an exact 
diagonalization using the Lanczos method and, there- 
fore, we can test the precision of the NCGSM/DMRG 
approach. Since details of the DMRG algorithm can be 
found in 6J], we only mention the basic ingredients of 
the DMRG calculation. 

It is important to consider in the reference space A, 
not only resonant states, but also some non-resonant 
shells from the B space, in order to generate all pos- 
sible many-body configurations in the warm-up phase, 
which these shells could generate 64] . In addition to the 
Osi/2 neutron-proton resonant states, the space A con- 
tained also the last d 5 / 2 and d 3 / 2 HO states. The choice 
of the shell to be included in the pole space A is arbi- 
trary and the results do not depend on which shell is 
considered. The mixture of positive and negative par- 
ity states assures that we will not miss any couplings in 
the warm-up phase. In the space B the shells are or- 
dered like: {is 1/2l ip 3 / 2 , ipi/2 , id 5 / 2 , id 3 / 2 , • • • }, where i 
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direct diagonalization 
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sweep =1 



max dim = 387 



sweep =2 



998 



sweep =3 



50 100 150 200 250 300 350 400 

N step 

FIG. 5. (Color online) Iterative process of the DMRG ap- 
proach for e — 1CP 7 and including only waves with £ — 2 
(s — pd). Both Lanczos and DMRG diagonalizations of the 
Hamiltonian are shown. 



denotes the scattering shell starting from 0. In the case 
of the ei-states i S (0,4) denotes the HO radial quantum 
number. 

Usually in the DMRG applications, it is decided from 
the beginning how many states of the density matrix that 
have the largest eigenvalue will be kept. This number is 
then kept fixed throughout the whole iterative process. 
In this work we will use the truncation scheme defined in 
(|2"2"1) . The smaller the e, the more vectors are kept. Usu- 
ally with e = I0~ 8 the exact result is reproduced. How- 
ever, even for a larger e (10 -7 or I0 -6 ) the agreement is 
very satisfactory, especially if the DMRG calculations are 
followed for more than two sweeps (see also Fig. [5]). This 
method is called the dynamical block selection approach 
because the number of states N p kept changes during the 
iterative process to satisfy the truncation limit (|2"2"|) . In 
Fig. [SJ we show the iterative DMRG process that in- 
cludes partial waves up to £ — 2. The number of steps 
starts from zero, denoting the first shell in the sweep- 
down phase. As discussed earlier, each step corresponds 
to the addition of a new s.p. shell. We notice a peri- 
odic pattern of pronounced oscillations in energy, that 
appear in the middle of each sweep, with their amplitude 
continuously decreasing with the addition of more shells 
until the convergence is reached. Finally, the energy has 
converged in the end of the third sweep. The truncation 
error for this calculations was e=10 -7 resulting in a max- 
imum number of vectors kept N p ^85 and a maximum 
dimension -D max ^1,200 of a matrix to be diagonalized. 
For a direct diagonalization in this basis, the maximum 
dimension would be: 387,998 in m-scheme and 96,883 in 
J-scheme. The exact result of the diagonalization is: 

£ cxact = -7.840 MeV, 
and the result obtained by the DMRG is 



-7.832 MeV. 



In the following, using the Vi ow k Argonne vis poten- 
tial, we will test the convergence of the DMRG method 
with respect to the truncation parameter e which con- 
trols how many vectors of the density matrix are kept 
at each step. The results are gathered in Figs. [5] and [7] 
Additionally, in Table HI we show the number of vectors 
that are kept for different values of the parameter e and 
the corresponding energy. Each point in Figs. [6] and [7] 
corresponds to the value at the end of the fourth sweep of 
the DMRG process. The error bar reflects the extremum 
values of the energy in the last sweep (see also Fig. [8]). 
We observe that with decreasing the value of e the error 
bar also decreases and almost vanishes for precisions bet- 
ter than 10 -5 , i.e. the DMRG result quickly converges 
to the exact result for e smaller that 10 -5 . 

The purpose of this exercise is to test the extrapolation 
properties of the NCGSM/DMRG calculations. In gen- 
eral, one would like to perform calculations with the best 
precision possible, i.e., e ~ 10 -8 , what is presently possi- 
ble only in the lightest systems. Extending the NCGSM 



9 



-2 

_ -3 
"5 -4 



OJ 

c 



Ss 
0) 

c 



> 

cu .4 



: (a) 

\ 


DMRG f— i 


\ 


fit 


: \ 


exact 




E = -8.041 ± 0.130 - 









10 2 10 3 10~ 4 10~ 5 10~ 6 10~ 7 10~ 8 10 s 



: (c) 1 


DMRG i— »-i 




fit — ; 




exact 




V E = -8.304 ± 0.081 ~_ 









10 2 10 3 10 4 10 5 10 6 10 7 10 8 10 9 



II 1 II 1 1 

: (b) 


- 

DMRG i — ▼ — i 




I fit 




exact 




E = -8.176 ± 0.101 - 









Iff 2 Iff 3 10" 4 10" 5 10" 6 10" 7 10" 8 Iff 9 



: (d) 


DMRG i — — i '- 




I fit — : 




exact 




\ E = -8.366 ± 0.030 









iff 2 iff 3 iff 4 io 5 10 6 wy 7 -\<y s io 9 



FIG. 6. (Color online) Convergence of 3 H g.s. energy as a function of the truncation error e of the density matrix (|22[) . The 
s.p. basis consisted of all states up to £ — 4. Panel (a) corresponds to a fit with four points, panel (b) with five, panel (c) with 
six and panel (d) with nine points. The fit function has the form: E = -E cx tr P + b ■ e c . For more details, see the description in 
the text. 
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FIG. 7. (Color online) Same as in Fig. [6] but omitting the low-precision e = 10 3 point. 



applications to somewhat heavier systems, one needs to 
choose a different strategy that makes the calculation 
feasible and, at the same time, gives a rather precise re- 
sult. The 3 H studies serve here as a testing-ground for 
this investigation. In Fig. |6]we show calculations with e 
ranging from 10~ 3 up to 10~ 8 . Panel (a) consists of four 



points that correspond to several values of e from 10~ 3 
up to 10 -4 . The result is almost 600 keV away from the 
exact one and the extrapolation is also of a poor qual- 
ity, lying at ~ 380 keV away. The function we use for 
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TABLE II. Vectors of the density matrix that are kept for 
different values of e. The corresponding matrix dimensions 
and many-body energies for 3 H are also shown. 



# vectors dimension energy (MeV) 



-2.5 

-3 

-3.5 

-4 
-4.5 

-5 
-5.5 

-6 
-6.5 

-7 
-7.5 



5.5 

450 



10 


" 3 5 


43 


-5.648 


5.0-1CT 4 8 


99 


-6.878 


2.5-10" 4 13 


109 


-7.396 


10 


" 4 22 


173 


-7.821 


5.0- 


10~ 5 28 


308 


-8.042 


10 


" 5 46 


600 


-8.287 


10 


" 6 73 


1075 


-8.357 


10 


" 7 96 


1909 


-8.381 


10 


" 8 117 


2575 


-8.388 




DMRG iteration scheme during the fourth sweep 
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FIG. 8. (Color online) Energy variation with respect to N s te P 
of the DMRG calculation for H and for several e values. The 
iteration is during the fourth sweep. 



fitting/extrapolating numerical data has the form: 

E = E ex t Tp + b ■ e c 

Adding one more point, the situation is improved, but 
still the exact result is off the error bars of the extrap- 
olation (panel (b)). In panel (c) the exact result is in 
between the error bars of the extrapolation and consid- 
ering all values up to 10 -8 we fall onto the exact result. 

The outcome is that carrying calculations with a rela- 
tively low precision, one can find the exact result within 
the error bars of the extrapolated value. In this specific 
example, the maximum dimension for e = 10~ 5 (panel 
(c)) is .D max = 600 (see Table HI)) . whereas the dimension 
in a direct diagonalization is 890,021 in m-scheme and 
123,835 in J-scheme. 

In Fig. [7]we make the same analysis but neglecting the 
energy point corresponding to e — 10 ~ 3 . This energy is a 
rather poor approximation of a final result, lying almost 3 
MeV away. By extrapolating, we retrieve the exact result 



even in the case (panel (a)), where energies that corre- 
spond to e up to 10 -4 were used in the fit. For this value 
of e, only N opt = 22 vectors of the density matrix were 
kept, resulting in a maximum dimension of the matrix 
to be diagonalized .D max = 173, which is much smaller 
than the dimension in a full diagonalization. Notice that 
the extrapolated value in panel (a) is not accompanied 
by an error estimate. This is because we are fitting here 
three points with a function that is characterized by three 
parameters and by definition the x 2 fit produces a zero 
error. 

In addition, in Fig. [3] we show the behavior of the 
energy during the fourth sweep in the DMRG process. 
Due to the small number of vectors kept, calculations 
that correspond to truncation parameters 10 -3 , 5.0T0 -4 
and to a lesser extent for 2.5T0 -4 , show variations with 
respect to the number of shells added (N step ). The error 
bars in Figs. [6] and [7] were calculated by the extremum 
values of these variations. For e values equal 10~ 4 and 
smaller, the energy is a flat curve with variations of the 
order of less than 1 keV. 

Finally, we have calculated the expectation value of 
the CoM operator (Eq. (HZ])). In the largest model 
space that we used, the expectation value of the H cm 
is approximately 7 keV. The fku energy in this case is 
18.5 MeV (b = 1.5 fm). Following the assertion of Ref. 
[6fj| . we conclude that in a sufficiently large model space 
the NCGSM wave function factorizes and the CoM wave 
function IV'Com) is a Gaussian with hcu — 18.5 MeV. We 
would like to mention that the cases of 4,5 He nuclei, re- 
garding the factorization of the wavefunction, is under 
investigation. 

The smooth convergence properties of the 
NCGSM /DMRG procedure both with the number 
of partial waves and the truncation error, as they were 
tested in 3 H and 3 He, will be used later to calculate 
somewhat heavier nuclei. One can perform several 
smaller scale NCGSM calculations, changing both the 
number of density matrix vectors and the number of 
partial waves, and extrapolate these data to retrieve the 
exact result within the error bars. 



C. He nucleus 

In this section we apply the NCGSM to 4 He nucleus, 
using the DMRG as the diagonalization technique. We 
will compare the NCGSM /DMRG results with FY calcu- 
lations, using the Argonne Ui§ interaction with a Vi ow k 
cutoff A = 1.9 fm -1 . The energies of the Osi/2 neutron 
and proton as they are calculated from the GHF process 
are -26.290 MeV and -24.453 MeV, respectively. As 
in the case of 3 H, the s and p scattering continua are 
taken along the real axis and are discretized with eigh- 
teen points. The contours extent up to 4 fm -1 . For the 
remaining d, /, g partial waves, we assume HO basis func- 
tions with b — 1.5 fm. We take 5ri, 3/ and ig states for 
protons and neutrons. For a truncation error e = 10~ 6 , 
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the maximum number of vectors, which are kept is N opt 
= 180, resulting in a Hamiltonian matrix of dimension 
-Dmax ~6,000. On the other hand, the dimension of the 
Hamiltonian matrix for a direct Lanczos diagonalization 
is computed to be 119,864,088 in m-scheme and 6,230,512 
in J-scheme. 

In Fig. [9] we show the iteration pattern for 4 He, which 
is similar to the one in 3 H case, but we observe that 
already in the middle of the third sweep the energy is 
converged. This can be attributed to the larger num- 
ber of vectors of the density matrix that are kept. The 
converged energy is: 



E 



NCGSM 



and the FY result [67j is: 



E 



FY 



-29.15 MeV, 



-29.19 MeV, 



Results of 3 H and 4 He using the Argonne vig are in a 
nice agreement with both FY and CC calculations with 
triples corrections [EH . We also performed calculations 
using the chiral N 3 LO interaction (A = 1.9 fm -1 ). The 
s.p. energies for the 0siy 2 neutron and proton poles are 
-24.333 MeV and -24.303 MeV, respectively. The first 
point in the number of steps corresponds to approxi- 
mately fortieth shell in the sweep-down phase. The re- 
sults are shown in Fig. 1 101 and the converged energy: 



E 



NCGSM 



= -27.48 MeV, 



is also in a good agreement with CC calculations with 
triples corrections [4JJ. 

The difference with the experimental binding energy is 
attributed to the missing three (or many) nucleon forces 
(3NFs) [H4zl- For a recent review of 3NFs we refer 
the reader to [7J]. In our work we neglected three-body 
processes, such as the two-pion exchan ge a mongst three 
nucleons, which appear at N 2 LO order [75| in the chiral 
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FIG. 9. (Color online) Comparison of the NCGSM with the 
FY result with the Vi ow k Argonne wis interaction. 
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FIG. 10. (Color online) Same as in Fig. |9]but with the Vi ow k 
N 3 LO interaction. 



Effective Field theory (EFT) framework [70. [77j|. They 
appear naturally in EFT together with a mid-range one- 
pion exchange and a short-range (contact) 3NF. Essen- 
tially, these parts are not included in our Hamiltonian. 
On the other, the fact that we employed a renormaliza- 
tion technique (i.e. Viowk) which integrates out high- 
momentum modes, will generate the appearance of 3NFs 
[M H M, M, M, [Zl- The latter are usually refercd to 
as "induced" 3NFs. 



D. 



'He nucleus 



5 He is a challenge for any many-body theory due to 
its unbound character. In particular, both the ground 
and first excited states are embedded in the continuum 
and they are particle unstable. They are many-body res- 
onances, which obey outgoing asymptotics. Because of 
these characteristics, the complex energy formulation of 
the shell model using the Berggren ensemble is very suit- 
able for its description. Contrary to the previous appli- 
cations, where we used a real-energy basis, in the case 
of 5 He we are employing a complex basis, which also in- 
cludes the 0p 3 / 2 resonant state (see also Fig. Q]and Fig. 
[2] for the position of this state in the complex /c-plane and 
its radial behavior). Overall we include the bound Osi/ 2 
neutron state with a s.p. energy of —23.290 MeV, the 
bound 0si/2 proton state with s.p. energy —23.999 MeV 
and the 0p 3 / 2 s.p. resonance with a real part of energy 
5Re = 1.193 MeV and a width V = 1267.12 keV. Its posi- 
tion in the complex-energy plane is: k = (0.277, —0.068) 
fm , The s.p. basis for protons and neutrons is pro- 
duced by a GHF calculation with an N 3 LO NN interac- 
tion, renormalized by the Vi ow k approach for a cut-off A 
= 1.9 fm" 1 . 

The p 3 /2 contour will be taken complex to satisfy the 
Berggren completeness relation (|11[) . whereas the s\/2 
and pi/2 contours may be chosen along the veal-k axis. 
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For states with I >1, we assume the HO basis functions 
(5d, 3/, 3g) as we described in the previous cases. 

In Fig. QT]we show the DMRG convergence pattern of 
the real part of the g.s. energy in 5 He. The calculation 
in Fig. [11] is presented starting from the fortieth shell in 
the sweep-down phase. The converged energy: 

K-Eshc = -26.31 MeV, 

lies at about 1 MeV above the experimental value and 
approximately 1.4 MeV below the energy found in CC 
singles and doubles approximation (CCSD) [4l|. The 
truncation error in this calculation is e = 10~ 6 and the 
maximum number of vectors we kept is N opt ~ 300. The 
corresponding dimension of the matrix is D max ~10 5 , 
whereas in the direct diagonalization one deals with a 
matrix of a dimension ^3xl0 9 . 

The imaginary part of the 5 He g.s. energy is shown in 
Fig. HH The converged value is: 

3£5 Hc = -0.2 MeV, 

i.e., rs He = 400 keV. For a comparison, the CCSD value 
for this width is Tccsd = 320 keV, whereas the experi- 
mental value is Fe xp = 648 keV. 

In the NCGSM, we are allowing all possible excitations 
of particles in the many-body wave function and include 
all possible nucleon-nucleon correlations generated by the 
NN interaction. Then a discrepancy between the results 
of NCGSM and CC calculations obtained for the same 
NN interaction is probably due to neglected correlations 
in the CCSD approximation. On the other hand, it is 
the chosen NN interaction that should be blamed for a 
deviation between the data and the NCGSM result. In 
this sense, the NCGSM result provides information that 
may help to improve the NN interaction and decide about 
the role of many-body interactions. 

In the DMRG calculation for 5 He, both the real and 
imaginary part of energy exhibit variations, which can 
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FIG. 11. (Color online) Same as in Fig. [I0]but for the real 
part of 5 He g.s. energy. 
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FIG. 12. (Color online) Same as in Fig. [11] but for the imag- 
inary part of He g.s. energy. 



be both negative and positive. The reason behind this 
behavior is that, in 5 He we use a complex basis resulting 
in a non-Hermitian complex symmetric matrix, thus the 
usual Ritz variational principle does not apply. On the 
other hand, in 3 H and 4 He, where the Berggren basis 
was real, the energy was gradually decreasing with new 
iterations. One should notice that even in the 5 He case, 
the real part of energy is almost always decreasing with a 
few notable exceptions during the first sweep, where the 
wave function is far from being converged. 

Having calculated 4 He and 5 He, one obtains the one- 
neutron separation energy 



Sin 



-1.17 MeV 



which, for a Vi ow k cut-off of A = 1.9 fm 1 , is in a rather 
good agreement with the data: S 1 ^ p = -890 keV. 



VI. CONCLUSIONS AND FUTURE 
PERSPECTIVES 

In this work, we applied the Berggren s.p. en- 
semble to perform ab initio NCGSM/DMRG calcula- 
tions for selected light nuclei, both well bound and 
unbound. We used a translational invariant Hamil- 
tonian and benchmarked our results against Faddeev 
and Faddeev- Yakubovsky calculations for 3 H and 4 He. 
We also investigated the extrapolation properties of the 
NCGSM/DMRG iterative procedure with the number of 
partial waves and the truncation error. We found that 
even if a relatively small number of vectors of the density 
matrix are kept, the NCGSM results can be extrapolated 
with high accuracy to the exact result. For heavier sys- 
tems, this methodology will be especially useful, since it 
will allow to mitigate the increase of the computational 
cost of the DMRG with increasing matrix dimensions. 

The NCGSM is a very natural choice to calculate un- 
bound nuclei, such as 5 He. For the description of 5 He, 
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we employed a complex-energy Berggren basis consist- 
ing of bound Osi/2 proton/neutron s.p. states, the 0p 3 / 2 
neutron s.p. resonance, and the associated real and com- 
plex non-resonant continua. We successfully reproduced 
the unbound character of this system from first princi- 
ples using the N 3 LO chiral potential as the NN inter- 
action. Further studies are necessary to understand a 
dependence of the energy and the width of 5 He, on the 
Viow k (SRG) renormalization scale parameter A (A) of 
the N 3 LO chiral interaction. 

The correct asymptotic behavior of the system and the 
coupling to the continuum plays an important role in the 
reaction theory. The GSM has been recently generalized 
for the study of reactions within a Coupled Channel (CC) 
GSM framework [Hi]. In this respect, the NCGSM can 
provide the realistic wave function for a target nucleus, 
which will include both many-body correlations and cou- 
pling to the continuum. 

In the near future, we plan to apply the NCGSM sup- 
plemented with the DMRG iterative procedure to calcu- 
late excited states of 4 He and 5 He, heavier weakly bound 
systems, such as 6 He, and very exotic systems in the hy- 
drogen isotopic chain. Understanding of the respective 
role of three-nucleon forces and a continuum coupling in 
light nuclei at the limits of nuclear stability will be an 



important challenge for future NCGSM studies. 
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